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Two conditional averages for the increment A{x,x + r) = T{x) — T{x + r) of the scalar 
field are estimated for DNS data: H{A) — <^ V|; A(a;, a; + r) | A(a;, a; + r) ) and G(A) = 
( I Vo; A(a;, a; + r)pl A(a;, a; + r)^ The probability density function -P(A) for A(a;,a; + r) is 
also calculated. If A, P(A), H{A), and G'(A) are expressed as 9, P{9), h{6), and g{e) in 
a dimensionless form, the data analysis has revealed interesting simple relationships between 
them: (1) g{9)P{9) ~ exp(— duh{u) / g[u)) derived by Ching and Kraichnan was numerically 
confirmed perfectly, (2) g{9)P{9) pass through an almost common point for any scale ranging 
from the dissipative to the inertial, to the energy containing scale, (3) by using fitting param- 
eters A and /3 as a function of r, h{9)/g(9) is satisfactorily fitted by (A//3) tanh(/3S), and then 
g{9)P{9) ~ (cosh/36')-'*'''' . These relations may be helpful to construct a turbulence model. 
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§1. Introduction 

Since Kraichnan's workjB the anomalous scaling of the 
passive scalar field has become a prototype.o£ a^^tudy of 
the intermittency problem of turbulence. EP'a-''u-''LP The 
scaling is evaluated by the scaling exponents of the nth 
order structure functions S'„(r) of the scalar field T{x)^ 
which are expressed as 



In eq.(h5) n is the molecular diffusion coefficient and x' 
is put equal to x + r. 



5„(r)-(A(a;,a; + r)"), 



(1.1) 



where A(a;, x + r) is the increment of the scalar field 
denoted by 



A{x,x') = T{x)-T{x'). 



(1.2) 



If the scalar filed is advected by the (5-correlated Gaus- 
sian velocity field, S'„(r) was shown to be governed by 
the following simple equation in the steady stateild^^B 

Dn{r) - Jn{r), (1.3) 
where Dn{r) is the inertial term expressed by 

d_ 

dr 



dr 



(1.4) 



where a(r) is the eddy diffusivity, and d, the dimension 
of the system. J„(r) is the dissipative term defined by 

J„(r) = nK{A{x, a;')"~^A(a;, x')), (1.5) 



where 



Aix,x') = iVl + \/l,)Aix,x'). 



(1.6) 
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In order to close the equation ( L3 ) , one has to make a 
certain approximation on Jnir). The ansatz, which was 
employed by Kraichnan,EP is that the conditional average 
of the field A{x, x') with a value of A(a;, x') being fixed 
is proportional to A{x,x') itself, i.e. 

H{A{x,x')) = {A{x,x')\A{x,x')) 

^ h{r)A{x,x'), (1.7) 

where fi{r) — J2{r)/ {2KS2{r)). Under this ansatz 
Jn{r) = nJ2{r)Sn{r)/{2S2{r)), which enables us to close 
the equation (1.3). AlthoughJthiSpansatz was examined 
in simulation and experiment, LP 'B-''lP'Ej) the definite con- 
firmation has not been reached yet. Hence, it is interest- 
ing to examine to what extent of A the linear relation 
(1.7) is satisfied. 

Another conditional average is conceivable. Making a 
partial integration in Jn{f) with respect to x, we have 

J„(r) = -n(n - l)K(A(a;. x + r)""^i;(a;, x + r)), 

(1.8) 

where 

I](a;, x + r) = \Vx{T{x) - T{x + r))|2. (1.9) 

Now we introduce another conditional average of 
I](a;, X + r) with a fixed value of A(a;, x -\- r) as 

G{A{x, x + r))^ {i:{x, X + r)\A{x, x + r)). (1.10) 

Although the linear ansatz has not been proved yet, 
it might have an origin in the dynamics of A(a;, a; + r). 
On the other hand, the conditional average ( 1.10| ) might 
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be closely related to the conservation of scalar disper- 
sion, since S(a;, x + r) represents the scalar diffusion 
rate. Which average is more appropriate will be judged 
through the examination of the result. This is the first 
aim of the present paper. 

The second motivation of examining those ajzHrages is 
the recent argument by Ching and Kraichnan,EiP which 
enables us to combine H and G to yield the probability 
distribution function (pdf) for A{x,x'). They showed 
that the pdf is related to the H/G in the homogeneous 
situation. We will confirm the relationship perfectly later 
in this paper. In this sense the conditional averages H 
and G are expected to shed light on the pdf of the incre- 
ment of the scalar field and the issue of intermittency of 
turbulence. 

Although the problem of the anomalous scaling should 
be ideally investigated in the randomly generated veloc- 
ity fields, the generation of such fields is very difficult. 
Thus, the above conditional averages are examined in 
the usual Navier-Stokes turbulence in the present paper. 

Finally we emphasize that the present paper is not 
directly devoted to the study of the recent issue of the 
anomalous scaling. One of the important features in the 
study is the introduction of the conditional average of 
spatial derivatives of the scalar field with a given value 
of its difference. Combination of the conditional average 
with the equation of motion and/or some kinematical 
constraints yields various check points to examine the 
pdf for the scalar difference, which would lead us to a 
more precise knowledge of the pdf. The point of view is 
not Imitted only to the case for the scalar field advected 
by the velocity field white in time. In this sense, we 
examine the statistics of tke scalar field advected by the 
Navier-Stokcs turbulencetJ' from the different point of 
angle with the help of the methods developed in the area 
of the anomalous scaling. 

The present paper is organized as follows. In the fol- 
lowing section the algorithm of the simulation as well as 
the statistical property of the velocity field is given. In 
§3 the pdf 's of the fields A and A are presented with em- 
phasis on their scale and Reynolds number dependence. 
To what extent of A the linear hypothesis of H is valid 
is examined in §4 and it is shown that it does not hold 
for large values of A, which correspond to the structure 
functions of orders larger than three. The deviation from 
the linearity is also discussed. In §5 the dependence of 
G on A is investigated; it is described by a simple fitting 
formula, which is valid up to structure functions of large 
orders. Section 6 is devoted to the pdf of A derived from 
the combination of H and G. The numerically obtained 
pdf agrees with the observed one perfectly. The remark- 
able findings about H/G and P{A)G are given and the 
results are analytically expressed using parameters. In §7 
the useful findings are summarized to help to construct 
a turbulence model. 

§2. Simulation 

The simulation was done with resolutions of 256^ (Run 
1 and 2) and 512^ (Run 3) mesh points for several values 
of molepilar viscosity using the NAL Numerical Wind 
Tunnel.Ej-' The velocity and scalar field are simulated at 



the same time in An cyclic box by making use of the 
Navier-Stokes equation and the advection equation. The 
Prandtl number is set equal to unity. The turbulence 
decays from the same initial spectrum for velocity and 
scalar fields 

^y|fc4exp(-2fc2). (2.1) 

The data were processed at t=10, when the energy dis- 
sipation rate reaches the largest value through the sim- 
ulation, but the scalar one has passed the largest value 
already. At that time the Taylor microscale Reynolds 
number R\ ranges from 89.5 for Run 1 to 120 for Run 3. 
The inertial region is identified from the scalar spectrum 
and the requirement of the relation T^Ur ^ r, where Tr 
and Ur are the scalar and velocity increment over r. The 
related data are summarized in Table 1. 

Table I. Statistical property of the present DNS. Run 1 and 2 
are carried on 256^ mesh points and Run 3 is on 512^ mesh 
points. Lp is the integral scale, A, the Taylor microscale, and rj, 
the Kolmogorov microscale. kmax is the maximum wavenumber 
after eliminating the aliasing error. 





Rx 




A/r, 


Lph 


inertial range 


Run 1 


89.5 


1.27 


18.6 


63.4 


7 <r/r) < 31 


Run 2 


120 


0.716 


21.5 


68.4 


12 < r/r; < 46 


Run 3 


120 


1.42 


21.8 


103 


14 < r/r; < 51 



Since the resolution for Run 1 and 3 are sufficient, the 
data associated with Run 1 and 3 are mainly discussed. 
The inertial region of Run 1 is located in 31 > r/77 >7 
for Run 1 and 51 > r/77 > 14 for Run 3, 77 being the 
Kolmogorov dissipation scale. The dissipation spectrum 
is the largest at about r/rj = 5 for both runs, so that 
the scale r/rj ^ 5, which will appear often later, is in the 
dissipative region. 

Since Run 2 and Run 3 have the same Reynolds num- 
ber with different resolution, the comparison of the result 
of Run 2 with that of Run 3 is interesting to know the 
effect of the difference of the resolution. 

§3. The pdf 's of the fields A^ and A^ 

In this section we will be concerned with the pdf 's for 
the scalar field increment Ar{x) = A{x, x + r) over two 
points separated by r and the Laplacian field increment 
Ar{x) = A{x,x + r). For simplicity we introduce the 
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normalized dimcnsionlcss quantities such as 

1 



(j)r{x) 



J2{r) 



Kr{x). 



(3.1) 



(3.2) 



The pdf's P{9r) and P{(j)r) of the fields 9r and (f>r are 
obtained for various values of r ranging from the dissi- 
pative scale to the upper inertial one. From now on, we 
suppress the subscript r without otherwise stated. The 
pdf's for Run 3 are depicted in Fig.l, although the same 
ones could be drawn for Run 1; r/rj — 20 in the iner- 
tial region in Fig. 1(a), and r/rj = 5 in the dissipative 
region in Fig. 1(b). The vertical quantities as well as the 
horizontal ones are normalized. The solid line denotes 
P{0), while the broken line does P{(f>)- The dotted line 
stands for the Gaussian distribution. In order to empha- 
size which part of the fluctuations contributes to the nth 
order structure functions, the arrows with P+ and P~ 
are inserted: P+ signifies a position where P{9)9^ is the 
largest for 6* > and P~ does a position where P{9)\9\'^ 
is the largest for < 0. 

The scale-dependence of pdf's of 9 and (j> is different. 
The frequency of large amplitudes of 9 becomes less as r 
increases, which implies that the degree of the intermit- 
tency measured from P{9) is weaken with the increase 
of the scale r. On the other hand, the pdf of (f> does not 
vary so much as compared with that of 9 even when the 
scale changes. This difference is consistent with the fact 
that 9r represents the fluctuations of scale r, while (j)r 
does the ones of dissipative scale. 

In order to show the Reynolds number dependence of 
both pdf's, we depicted P{9) at r/r] = 20 for Run 1 
and Run 3 in Fig.2(a) and P{(f>) in Fig.2(b). The sohd 
line is for R\ — 89.5 of Run 1, and the broken line is 
for R\ = 120 of Run 3, together with the dotted line 
for Gaussian. The large components of 9 and (j) become 
more frequent as Rx increases. 

When one compares P{<j)) to P{9) of the same size r, 
for instance by inspection of Fig. 1(a), one notices that (f) 
is more frequent than 9 for large amplitudes and less fre- 
quent for intermediate ones for any r. As r increases, the 
discrepancy between (p and 9 becomes large, although 
figures associated with the other scales are not shown 
here. We suspected that the difference between both 
pdf's would be responsible for the deviation of H from 
the linearity ansatz, which will be discussed in the follow- 
ing section, but we could not have found a clear reasoning 
for it so far. 

The Run 2 and 3 have the same Reynolds number 
R\ — 120, but are simulated on the different number of 
mesh points; the former is on 256"^ and the latter, on 
512^. Thus, the small scale resolution is less for Run 
2 than for Run 3. P(0)'s for both runs are almost the 
same in the inertial region, but P(0)'s are considerably 
different from each other even in the inertial region. The 
curve of H{A) for Run 2 is, however, found to be quite 
the same in the inertial range as Run 3. Hence the lin- 
earity ansatz (1.7) may be more robust in the inertial 
region than the pdf's. 




PI p' p'l p', PI 



Fig. 1. The pdf's of 6 and for Run 3: (a) r/r; = 20 in the 
inertial region, and (b) r/rj = 5 in the dissipative region. The 
vertical quantities as well as the horizontal ones are normalized. 
The solid line denotes P{Q), while the broken line does P(<t>). 
The dotted line stands for the Gaussian distribution. The arrow 
with Pn signifies a peak position where P{9)9" is the largest for 
6 > 0, while P^ is a peak position where P{6)\6\" is the largest 
for e < 0. 



§4. Examination of the linear ansatz 

The conditional average H is expressed in the dimen- 
sionless form as 



h{9r) = H{Ar 



J2{r) ■ 



(4.1) 



Note that h{9) is positively defined for 6 > because 
J2(r') is negative. Under this definition the linear ansatz 
(1.7) becomes 



h{9) = 9. 



(4.2) 



We investigated h{9) for various values oir/rj for every 
run. The observed results are summarized. The linear 



relation (4.2) is approximately valid, but not in detail. 
Figure 3 is the curve of h{9) against 9 ior r/rj — 20 (the 
solid line) and r/rj = 30 (the broken line) for Run 3. 
Note that the same type of curve holds for any inertial 
range separation for other runs. The dotted line stands 
for h{9) = 9. In the inertial region h{9) is slightly larger 
than 9 for small values of 9, while for larger values of 9 
h{9) is considerably smaller than the linear relation as 
seen in Fig. 3. Within the present simulations we did not 
see the Reynolds number dependence of h{9), although 
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R,=S9.5 - 

R,=no 
Gaussian 




(b) / 


\ /?.=89.5 




\ i?,= 120 - 




\ Gaussian 







Fig. 2. The Reynolds number dependence of pdf's of 6 and ip at 
r/rj = 20 in the inertial region: (a) P{6) and (b) P{<j>)- The 
vertical quantities as well as the horizontal ones are normalized. 
The solid line is for R\ = 89.5 in Run 1, and the broken line is 
for R\ = 120 in Run 3. The dotted line stands for the Gaussian 
distribution. 
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r/7]=20 — 
r/T)-30 — - 
Ansatz 
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Fig. 3. The plot of h{e) against 9 at r/r; = 20 (solid) and 30 
(broken) for Run 3. The data in the inertial region drop on 
nearly the same line. The dotted line stands for h{9) = 9. The 
arrows with Pn have the same meaning as before. In the inset 
h{9) is plotted for the dissipative range separation r/rj = 5. 



Ching, L'vov and Procacciai^ observed by processing the 
experimental data that the hnearity would become more 
appropriate as the Reynolds number increases. 

In Fig. 3 the arrows with P„ are inserted for conve- 
nience, so that the nth order structure functions are pre- 
dominantly contributed by the components around P„. 
The statistical reliability is assured up to Pq or P-j at 
best. Note that h{9) is closer to the linear ansatz about 



P2, which is confirmed for any separation in the inertial 
region. On the other hand, the linear relation is not good 
around P„ with n > 3. This suggests that the values of 
Sn with n > 3 do not rely on the linear hypothesis in the 
present simulation. 

The function h{6) must satisfy the conditions 



and 



{h{e) - d)9P{9)d9 ^ 0, 



{h{9)e -i)p{e)de = o. 



(4.3a) 



(4.3b) 



which are derived as follows. The definition of H yields 

J2 
2k 



J AH{A)P{A)dA = (AA) = -(E) 



whose dimensionless form reduces to / 9h{9)P{9)d9 = 1. 
On the other hand, 9 is normalized so that / 9^P{9)d9 — 
1. Combi ning these two equations, we are led to the 
condition ( 4.3a ) a nd (4.31: ). 

The condition (4.3a) is important in the following 
sense. Since P is a functional of h{9) and positive defi- 
nite, the function h{9) — 9 must change its sign at least 
once if there exists any deviation of h{9) from the lin- 
ear ansatz. Function 6P{9) decays rapidly at large 9 so 
that small deviation of h{6) from the ansatz at small 9 
implies large deviation of h{9) at large 9. Equation (4.3) 
also can be used as an examination of the linear ansatz 
in that smaller deviation of h{d) — 9 from at small 9 
where the sample size is usually enough leads to smaller 
deviation at large 9 where the sample size is poor. This 
is consistent with observation in Fig. 3. 

In the dissipative range with r /"q = 5, h{9) is given in 
the inset of Fig. 3. The function h{9) is slightly smaller 
than 9 for small values of 9, and larger than one for larger 
values of 9 in contrast to the inertial range separations. 
The crossing of h{9) with the linear relation 9 occurs 
around P3 for r/rj = 5. This tendency is observed also 
in Run 1. 

In order to know how much h{9) is deviated from 9 we 
evaluate a variance 



{h{9) - 9fP(9)d9 



(4.4) 



for run 3. At r/?7 = 5 in the dissipative region a — 0.008. 
As r increases, it decreases and reaches zero at r/rj — 10, 
where the linear ansatz is good. As r increases further 
in the inertial range, a grows; a — 0.01 at r /r] = 30. 

The observation that h{9) deviates from the linear re- 
lation for n > 3 in the inertial region is understood in 
the following way. The ansatz (1.7) is the projection of 
the field K{x^x') on the field A{x,x'). The coefficient 
/i(r) is obtained by multiplying eq.(1.7) by A{x, x') and 
taking the average over the distribution of A{x, x'): 



h{r)S2{r) = {A{x,x')H{A{x,x'))) 



{A{x,x')k{x,x')) 



■Mr) 
2k ' 



(4.5) 



from which we have fi{r) ~ J2{r) /{2KS2{r)). The cal- 
culation of J2{r) is the same as that of the projection 
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coefficient. Since J2{r) is directly related to S2{r), it is 
naturally expected that the linear ansatz including the 
definition of fi{r) is well satisfied at n = 2, i.e. around 
P2. Let us go to larger values of n. In the calcula- 
tion of J„(r) the field A{x,x') must be projected on the 
space where the quantity A(a;, £c')"~^ is dominant. Such 
a space is considerably localized, so that the projection 
coefficient /i(r) is not suitable. /i(r) must be replaced 

by 



Mr) 
2KS2{r)' 



(4.6) 



where J2(?') and 5*2 (r) are computed over such an en- 
hanced space. Since the equation for the scalar field is 
linear, one may expect that the ratio of J2{r) to S2{r) is 
constant everywhere. This would be correct if the scalar 
field did not move in space and the cascade occurred lo- 
cally in space. While the strong fiuctuations of scale r 
cascade down to the dissipation scale where the field A is 
mainly contributed, they will diffuse to the less intense 
space through the convection term due to the velocity 
field and the viscous term due to the molecular viscosity. 

Although J2{r) is larger than J2{r), a ratio of J2{r) 
to J2{r) is less than that of S2{r) to 5*2 (r). Hence 
is expected to be less than fi{r), implying that h{6) is 
lower than 9. On the other hand, weak fiuctuations will 
be affected in the opposite way, which may be responsible 
for the positive deviation from the linearity. 







/ 






\ -6h(B) 




\\ g(6) 







Fig. 4. The conditional average of the spatial transfer of the en- 
ergy flow, i.e. t{e) = g(e) - eh(e) for r/r; = 10 in Run 3. The 
solid line stands for t((9), the broken one does for —6h{6), and 
the dotted one does for g(6). A positive value of t{9) corresponds 
to the inflow and a negative one does to the outflow. 



The above interpretation can be supported more di- 
rectly by investigating the conditional average of the fol- 
lowing identity with a fixed value of A{x, x + r): 

A{x, X + r) A(a;, x + r) ^ 

Vx ■ {A{x,x + r)\/xA{x,x + r)) - |Va;A(a;, a; + r)p. 
Using the notation introduced previously, we have 
(Va, • {A{x,x + r)VxA{x,x + r))\A{x,x + r)) 

= A{x, X + r)H{A{x, x + r)) + G{A{x, x + r)), (4.7) 



where the left hand side is the conditional average of the 
spatial transfer rate with fixed A{x,x'). Since the two 
conditional averages on the right hand side are numer- 
ically evaluated, the dimensionless conditional transfer 
rate t{d) (= g{9) — 9h{9)) is calculated as given in Fig. 4. 
The region with positive values of t{9) is characterized by 
the inflow, while that with negative values is by the out- 
flow. About P2, t{9) ^ as expected. For large values of 
A{x,x') the outflow is dominant. Hence /i(r) < /i(r) is 
anticipated in such a region. For small values of A{x, x') 
the inflow is dominant, so that the opposite tendency is 
expected. 

§5. Examination of the second conditional aver- 
age 

We studied G{Ar), the conditional average of S(a;, x + 
r) with A{x,x + r) fixed. For simplicity we express G 
in the dimensionless form 

2k 



Mr) 



G(A.). 



(5.1) 



The analytic expression g{9), if possible, is very useful 
to discuss the pdf P{9), since it is directly related to the 
pdf as shown in the following section. Hence we try to 
investigate the funtional form of g{9) in detail. 

First we plot g{9) against t^. The reason for it is that 
there must be a certain conservation relation between 
the scalr dispersion 9^ and the scalr diffution rate g{9) 
for large values of 9. Figure 5 is such a plot for \9\ < 5. 
According to it, g{9) may be fitted by a quadratic form 



m ^ goir) + gi{r)9' 



(5.2) 



for large values of 9. For r/rj = 5 in the dissipative region 
the curve is convex, so that ( |5.2| ) is not valid; g{9) ^ 9" 
with a < 2. In the inertial region r/rj = 1 5 a nd 20 
we can see the considerable range fitted by ( |5.2| ); the 
dotted straight line drawn on the curve for r/r] — 15 is 
expressed by 1.7-1- 0.126*^. The linearity is remarkable in 
the range between P3 and P5 (which are calculated for 
r/rj — 15), where the linear ansatz (1.7) fails, as seen in 
Fig. 3. Such a fitting ra nge becomes short as r increases. 
By fitting the formula (5^) to the curve oi g{9) of Fig. 5 
in the range 7 < 9^ < 13, we estimated go{r) and gi{r) 
as a function of r; they are depicted in Fig. 6; the sohd 
line is for go and the broken one is for lOgi. Both are 
slightly decreasing with r. Note that gi{r) is one order 
of magnitude less than go{r). 

In the region with small values of 9 g{9) deviates from 
a fitting formula (5.2). To clarify the point, we depicted 
g{9) against 9 in Fig. 7 in the range \9\ < 2. This figure 
indicates that the curve is described by a quadratic form 
in the vicinity of ~ and then, by the linear form for 
larger values of 9 around 1. The linear region may be 
a transition region between the two different quadratic 
forms. 

Before concluding this section we want to add that the 
function of g{9) must satisfy the relation 



{g{9) - l)P{9)d9 = j {g{9) - 9^)P{9)d9 = 0, (5.3) 
which is derived in a similar way to (4.3a). From the 
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r/lj - 5 - 
//IJ- 10 - 
r/lj- 15 
r/lJ - 20 
r/lj - 25 
r/lJ -■ 30 - 
r/lJ -■ 35 
/•/IJ -■ 40 



Fig, 



5. The plot of 3(0) for various r for Run 3. 
t)2 



It is plotted 

against 6^. The linearity is quite good for r/rj = 15 as seen 
from the dotted straight line 1.7 + 0.12(^^ drawn on the curve for 
r/rj = 15. 



lOg, 



r/T] 



Fig. 6. The plot of go{r) and gi{r), which are estimated from 
Fig. 5 in the range 7 < 6^ < 13, as a function of r/r; in the inertial 
region for Run 3. The solid line is goir), while the broken line is 
Wgi(r). 
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^ r/l)~10 




• r/l7~15 




/ r/l) ~ 20 




/ r/l) ~ 25 ■ - - - - 




r/lJ - 30 




r/T] ~ 35 




/•/I) ~ 40 



Fig. 7. The plot of g{6) for various r for Run 3. 
against 6 in the range of small values of 8. 



It is plotted 



definition of G, its dimensionless form satisfies 
g{0)P{0)d0 = 1. 



Combining it with / P{e)d9 = J e^P{e)de = 1, we have 
the condition (|5.3D. It is guaranteed in Fig. 5 and 7. 



§6. Derivation of the PDF of 9 in terms of h{e) 

and g{e) 

Ching and KraichnanS^ argued that P(6) in the ho- 
mogeneous system is expressed in terms of h{d) and g{9) 
through 



P(%(0) = Cexp 



du q(u) 



where C is a normahzing constant and 

h{u) 



q{u) = 



(6.1) 



(6.2) 



Since the derivation of the relation is straightforward 
and short, we repeat their derivation for convenience. 
The extension of (4.7) to any n yields 

(Va; • (A(a;,a; + r)"-iVa;A(a;,a; + r))|A(a;,a; + r)) 

= A"-i(a;, X + r)H{A{x, x + r)) 

+ {n - l)A"-2(a;^ X + r)G{A{x, x + r)) (6.3) 

If the system is homogeneous in space, the average on 
the left hand side will be replaced by the spatial aver- 
age, so that the left hand side vanishes. Then in the 
dimensionless form this equation becomes 



dee''-'^h{e)p{e) = {n - 1) / dee'^-^g{e)p{e) 



9e^'-'^g{0)P{e). (6.4) 

The above relation must hold for any n, so that 

mm - -^g{e)P{e), (6.5) 

whose solution is ( ^.l| ) . Hence the relation is expected to 
be exact for a spatially homogeneous system. Note that 
the similaiL^xpression for the pdf was proposed already 
by YakhotllJ) in Benard convection. 

First, we present the observed result of P{0)g{9) for 
Run 1 in Fig. 8(a) and for Run 3 in Fig. 8(b) for various 
values of r ranging from r/rj — b'm the dissipative region 
to r/rj = 100 in the energy containing region. The Gaus- 
sian pdf is also added. (1) The product form P{9)g{9) 
looks simpler than P{9) itself in Fig.l. The fitting for- 
mula for P{9)g{9) will be given later. (2) All curves pass 
through an almost common point. According to Fig. 8(a) 
and (b) the point is read as {9 = 1.5, P{9)g{9) = 0.13). 

When we magnified the crossing region of the curves, 
however, we noticed that the crossing does not occur at a 
single point. All the curves for the inertial region ranging 
from r/rj — 15 to 50 (Run 3) cross at 6* = 1.31. On the 
other hand, the curves for the energy containing region 
and the dissipation region meet those for the inertial re- 
gion at shghtly larger values of 9, i.e. about (1.4 ~ 1.5). 
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Second, we pay attention to the function q{9) = 
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Fig. 8. The product of P{Q) and g{d) for various values of r/r?. 
The Gaussian pdf is also inserted, (a) Run 1 and (b) Run 3. 



Fig. 9. The function q(d), which is numerically evaluated with 
the aid of the observed h(9) and g{9) for Run 3, is plotted against 
9 for various values of r/ri. 



Si' 









X rhsof(6.1) 


(a) / 









h{9)/g{9). Substituting the measured values of h{d) and 
g{6) into q{0), we numerically obtained q{9) as a func- 
tion of 9; the function q{9) for various values of r/rj for 
Run 3 is given in Fig. 9. It increases linearly from the 
origin. The slope decreases with increasing r, but the 
decrease is not appreciable in the inertial region. For the 
separation in the dissipative region, q{9) saturates and 
then, decreases very slowly. As r increases, it does not 
saturate enough. The curves lose the symmetry a little 
when \9\ goes beyond 4. When r increases, the consid- 
erable fluctuations are observed in the region with large 
values of 9. Since q{9) looks simpler than h{9) and g{9) 
themselves, the fitting formula for q{9) is more suitable, 
as developed later. 

Since the simulated system is co nsidered to be homo- 
geneous in space, the relation (3J) is expected to hold. 
Su bstit uting the calculated q{u) into the right hand side 
of ( xl ) and carrying out the numerical integration over 
u, we evaluated the right hand side of (6.1). Then the 
result was compared with the observed P{9) multiplied 
by g{9). The comparison was shown for r/rj — 5 and 
r/rj = 20 of Run 3 in Fig.lO(a) and 10(b). The agree- 
ment is perfect. 

Before going to discuss the analytical form for q{9) and 
g{9)P{9), we would like to consider which part of fluc- 
tuations contributes to q{9), and g{9)P{9) consequently. 
Let us start with small values of 9, where h{9) —9 and 
g{9) = const = cq. Hence, q(u) — u/cq , leading to the 
Gaussian distribution according to (3.1) and the curva- 
ture is inversely proportional to cq. Since cq decreases 





p{e)g(e) — 




\ rhsof(6.1) 






/ 





Fig. 10. The comparison of the numerically integrated right hand 
side of eq.(6.1) with the observed pdf multiplied by g{6) for 
(a)r/77 = 5 and (b) r/r] = 20 for Run 3. 



as r decreases, as seen from Fig. 7, the curvature around 
= is a decreasing function of r, in agreement with 
Fig.8. 

The deviation from Gaussian begins to occur as 9 
grows. As 9 increases, h{9) is still close to 9, but 17(6') 
is not constant. Figure 7 suggests that g{9) is approx- 
imated by g(9) = cq + Ci9 in that region, where cq in- 
creases and ci decreases with increasing r. Then 

g{9)P{9)^Ce.J-^+'^log'-±^). (6.6) 
V ci c| co/ci / 

The existence of ci makes the pdf deviate from Gaussian. 
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If the first term is dominant, the pdf is of the exponential 
type: the exponent 1/ci increases with r, in consistent 
with the observation of Fig. 8. 

To conclude this section, we would like to propose an 
analytic expression of q{0). The analytical expression, if 
possible, is helpful for constructing a turbulence model. 
The combined function q{9) is simpler than h{9) and g{9) 
themselves. It may be approximately parametrized as 

q{e) = ^ts.nhm, (6.7) 

P 

where A and /3 depend on r. The most optimal fitting 
values for them are depicted in Fig. 11; the fitting was 
done in the symmetric region \d\ < 4. The crosses stand 
for P, while the pluses do for A; A and /3 are decreasing 
with r, while A/P is increasing. Although the best fit 
for the entire region of is given in that figure, the local 
fitting is also possible. (For instance, the separate regions 
9 <1 and 9 > 1 would be fitted more appropriately.) 
Substituting (6/7) into (|6.l|) yields 



P{9)g{9) = C(cosh p9)-^^>^' , (6.8) 

where C /3/S(A/2/32, 1/2), with B being a beta func- 
tion. For smaU 9, P{9)g{9) - e--^^'/^, i.e. Gaussian. 
The curvature decreases as A when r increases, in agree- 
ment with Fig. 8. For large 9, it becomes e~^^/^ , i.e. the 
pdf is of the exponential type, whose coefficient is A/ f}, 
which is an increasing function of r, also consistent with 
Fig.8. 

In order to show how the fitting works, we would like 
to give g{9)P{9) for r/77 = 10 as an example in Fig. 12, 
where A — 1.44 and (3 — 1.60 for r/rj — 10 are employed 
from Fig. 11. The solid line stands for the observed value 
of g{9)P{9) in Fig. 8(b) and the broken line does for the 
right hand side of (3.1) obtained after integrating q{u) 
over u. The dotted line is given by ( ^.^ ) with selected A 
and (3. The agreement is very satisfactory. 



CO. 



A 



5 10 15 20 25 30 35 40 45 50 



r/7] 



Fig. 11. The fitting parameters A and /3 as a function of r for 
Run 3; the crosses stand for /3, while the pluses do for A. 
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Fig. 12. The comparison of ( p.^ ) to the observed function of 
g{8)P{8) at r/rj = 10. The solid line stands for the observed 
value o f a(9 )P(9) and the broken line does for the right hand 
side of (6.1) obt aine d after integrating fi'{'ii) over it. The dotted 
line is given by (Oh with A = 1.44 and /3 = 1.60. 



§7. Discussion 



When the equation (1.3) is expressed through the pdf, 
it takes the form 



->^-g^{P{A,r)H{A,r)) 
-n^{P{A,r)G{A,r)). 



(7.1) 



If the velocity field is not 5-correlated in time, the inertial 
term will be more complicated. In any case, however, 
P(A), H{A) and G(A) must satisfy similar equations. In 
this sense the information on H and G extracted from a 
study of the simulation is useful to construct a model for 
turbulence. For this purpose we summarize the useful 
points. 

1. The combined form q(9) is a more appropriate func- 
tion for modeling, rather than h{9) and g(9) them- 
selves. 

2. Another combined form P{9)g{9) is more appro- 
priate than the pdf P{9) itself. P{9)g(9) is the 
weighted conditional average of the scalar diffuison 
rate g{9). 

3. The curves P{9)g{9) with varying r drop on an al- 
most common point. 

Let us consider the meaning of g{9). Take a specific 
scale r. Even if 9r is fixed, there are many structures, 
so that the scalar diffuison rate \Wx9{x,x + r)p vary 
from structure to structure. After taking the average 
over structures, we have the conditionally averaged dis- 
sipation rate g{9r)- 

The gr{9r)Pr{Sr) is the average diffusion rate at scale 
r weighted by the frequency Pr{9r) of 9r- Here we insert 
the subscript r for the functions g and P to show their 
scale dependence explicitly. The present simulation in- 
dicates that gr{9r){= gr{9r)Pr{9r)) takes the same value 
0.13 at 9r = 1.5 irrespective of r. For 9r < 1-5 gr{9r) in- 
creases with increasing r, while it decreases for 9r > 1.5. 
This is surprising, we think. Why? To show it, we con- 
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sider two specific scales ri and r2 such that ri < r2- 
Then we can have a situation that {d) < Qr^ {0) for 
9 < 9c, where 9c is a certain value, while gri{d) > SrsC^) 
for 9 > 9c. However, the present DNS suggests that this 
inequality occurs at the same 9c = 1-3 for any inertial 
separation r, ranging from the dissipative to energy con- 
taining scale. 

Finally we speculate the pdf in turbulence for large 
Reynolds numbers based on the results, although the 
present DNS decaying turbulence has only a limited in- 
terval of the inertial range. In this paper we observed 
the region where g{9) oc 9^ for large values of 9. If such 
a term is dominant, and h{9) is 9, the function q{9) is 
proportional to 9~^, yielding the omver law distribution 
as suggested by Sinai and Yakhot,E3' and Ching.ES Fig- 
ure 5 indicates that the dominance of the square term 
over the constant term in g{9) becomes more enhanced 
as r decreases, so that the pdf of the smallest scale is 
expected to become more power-law like. However, in 
this simulation the pdf for the smallest scale is still of 
exponential type. The reason for it is that h{9) ^ 9, but 
g{9) is weaker than 9^, as seen in Fig. 5, so that q{9) at 
r/r] = 5 decreases slightly for large 9 (not as 1/9), as 
seen in Fig. 9. In order to know the behavior of the pdf 
for large amplitude, therefore, the information of q{9) is 
more indispensable than that of individual conditional 
average h{9) and g{9). 
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